SYSTEM IDENTIFICATION FROM 
NONUNIFORMLY SAMPLED DATA 


A thesis submitted 

in partial fulfilment of the requirements 
for the degree of 


MASTER OF TECHNOLOGY 


By 

ABHAY C RANADE 


* «» ! 


to the 

JEPARTMENT OF ELECTRICAL ENGINEERING 
SDIAN INSTITUTE OF TECHNOLOGY KANPUR 


FEBRUARY 1990. 



CERTIFICATE 




:^ DrMr Sv 


.o 





1 f . 88 
I'btnitti * 




It is certified that this work entitled SYSTEM IDENTIFICATION FROM 
NQNUNIFORMLY SAMPLED DATA by ABHAY C RANADE has been completed under my 
supervision and that this work has not been submitted elsewhere for the award of 
a degree 


IA I r Z 111° 

Dr P Sircar / ' 

Assistant Professor 
Dept of Electrical Engineering 
Indian Institute of Technology, 
KANPUR - 2D8016, INDIA 



- 9 APR 1990 

CENTRAL LIBRARY 

I ! ~ KANPUR 


' 5 


Acc No 


Al()7flSo 

“ • I ■ * 4 


£_£' l^qo-M- fV^N- S^2 




f ♦ 



ACKNOWLEDGEMENT 


It is my pleasant duty to express my sincere gratitude to Dr P Sircar 
his encouragement, valuable comments and constructive suggestions 1 would lik 
thank him for his many ideas which helped me to broaden my perspective on si r 
processing 1 express my appreciation to him for guiding my efforts and sbar 
my enthusiasm in this challenging project 

No words are adeciua te to express by indebtedness to my mother for all 
pains and sufferings she has undergone to bring me up to this stage 


i 

i 


l 


> 


I 






ABSTRACT 


In this thesis, a technique has been presented for system identification 
from a nonumformly sampled data In this regard, the idea of Orthogonal 
Polynomial Approximation is used for signal reconstruction Acourate parameter 
estimation is achieved by using the statistical properties of the colored noise in 
the reconstructed signal The estimation procedure presented here makes use of 
the subspace separation approach Singular Value Decomposition (SVD) has been 
used, in this regard Furthermore , the ever present problem of model order 
selection has been addressed adequately 

The technique presented here, was also applied to the special case of 
uniform sampling Performance comparison, with a conventional technique, for the 
case of uniformly sampled data, proved the superiority of the technique given 
here 
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CHAPTER 1 
NTROOUCTTON 


Accurate estimation of parameters of a signal embedded in noise is a 
problem faced in many signal processing applications Particularly important is the 
case when the parameters to be estimated are time-invariant These parameters are 
estimated from a given set of measurements of a signal These measurements can 
be performed in two ways, viz , by uniform sampling or by nonuniform sampling In 
this thesis we consider the case of nonuniform sampling 

Nonuniform sampling can be chosen over uniform sampling for various 
reasons Efficiency of sampling and noise reduction are two such considerations 
favouring nonuniform sampling 

An efficient sampling strategy would be one which has a capability of 
adjustable sampling rate The sampling rate is adjusted, according to the rate of 
change of signal Thus, more number of samples are taken in the region where the 
signal changes rapidly, and reduced number of samples in the region where the 
signal varies slowly 

Yet another sampling strategy will be to use nonuniform sampling to 
increase the SNR of the sampled points This can be achieved for a signal 
embedded in zero mean noise, by sampling more number of points when the signal is 
away from zero crossing, and by sampling less number of points when the signal is 
near the zero crossing The reasoning behind the strategy should be apparent 

In both of the sampling schemes, it is evident that the sampling points will 
be distributed arbitrarily, although it will be necessary that the sampling theorem 
be satisfied locally Cil in order to aohieve complete reconstruction of the signal 

Modern signal processing methods use parametric modeling to achieve high 
resolution and accuracy The idea is to choose a particular model depending on a 
priori information about the signal Use of incorrect model, however, deteriorates 
the overall performance 

Modeling of an ARMA prooess by an AR process is a typical example, where 
use of incorrect model reflects in certain fundamental limitations like high model 
order reauirement and increased computational burden Another example is of 



parameter estimation of a signal consisting of damped sinusoids Use of adaptive 
spectral estimation techniques in this case has been suggested In such an attempt, 
however, the model itself becomes Ume-varying Thus, it becomes necessary to 
adapt the model parameters to the nonstationanty m the given signal Needless to 
say, such methods are restricted to only those applications where a large record 
of data is available, and also the damping factors are not too large 


In this thesis, we present a time-invariant model to characterise the signal 
consisting of complex exponentials embedded in noise Thus, our requirement is to 
estimate a constant set of parameters accurately, instead of fast and reliable 
tracking of time-varying parameters Note that the use of small data records for 
parameter estimation is our greatest constraint, and under this constraint we will 
fully utilize the nonuniform sampling scheme, and achieve efficient sampling 
together with high accuracy in estimation However, throughout this thesis we will 
maintain arbitrariness of the distribution of sampled points, except that the 
condition for complete reconstruction of signal is satisfied 

Sircar and Sarkar Cil have worked on system identification from a 
nonuniformly sampled signal It has been shown that the estimated poles are 
accurate only for high SNR, but they deteriorate substantially as the SNR is 
reduced, a factor which is attributed to the numerical ill-conditioning of the 
problem itself Neverthless, the fact that their model dose'nt take into account the 
embedded noise should not be overlooked Thus, by utilizing information about the 
noise while estimating the signal poles, it oan be expected that substantial 
improvement in accuracy can be obtained 


In this thesis, we make use of the so called Orthogonal Polynomial 
Approximation, as used by Sircar and Sarkar lil, for reconstruction of 
nonuniformly sampled signal In the present work, we will also characterise the 
noise present in this approximating (reconstructed) signal and use that information 
to get accurate estimates of the signal parameters 

It was observed m the process of simulation, that the noise present in 
the approximating signal was colored, even when the noise embedded m the signal 
was white 



A colored noise sequence is characterised by the fact that it is 
correlated, unlike white noise sequence, which is uncorrelated Thus, colored noise 
is a nondetermimstio, bandhmited signal 

A major contribution of this thesis is that a technique has been 
suggested which can estimate the signal parameters even in the presence of 
colored noise at the first place In this regard, we assume that the so called 
Normalised Covariance Function 'of the noise, embedded m the nonuniformly sampled 
signal, is known a prion For the purpose of simulation we assumed the noise to be 
white, Gaussian and Wide Sense Stationary (MSS) 

It is interesting to note that the colored noise in the approximating signal 
originates because of the polynomial approximation procedure Note that the white 
noise assumption is not stringent, since, once the Normalised Covariance Function 
is known, both the cases, viz , white or colored, become exactly same 

A variety of applications which can use the above mentioned methodology, 
are given below 

In Electro Magnetic Pulse (EH 3 ) applications [23, one of the primary 
requirements is to obtain the transient response of a test object In this scheme, 
an Electro Magnetic Pulse is impinged on the test object The output of the object 
is sensed by a probe, having a known transfer function, and is digitised These 
samples are then used to identify the transfer function of the test object A 
typical EMP has a short rise time and rather slow decay rate Thus, nonuniform 
sampling will certainly be bemf icial in this case 

Another application could be in high resolution Direction Of Arrival 
(DOA) estimation C33 DOA estimation is important in many sensor systems such as 
radar, sonar, electronic surveillance and seismic exploration The transmission 
medium is assumed to be isotropic and nondispersive so that the radiation 
propogates in straight lines, and the sources are assumed to be in the far-field of 
the array Me know that if and only if the transmission medium is vaccum, the 
signal propogates undamped But m applications like radar, sonar and seismic 
exploration, the signal could be sufficiently damped, thereby making it 
nonstationary [43 Algorithms like MUSIC (Multiple Signal Classification) [53 and 
ESPRIT (Estimation of Signal Parameters via Rotational Invariance Technique) C63 



have been suggested for solving the above mentioned problem MUSIC assumes that 
the signal is stationary, which may not be the case in true sense Also, both MUSIC 
and ESPRIT use uniform samples and assume that the embedded noise is white The 
later assumption is not true in most of the practical cases 

A very important application is m Linear Predictive Coding (LPC) I?D for 
bandwidth reduction of speech signals The idea is to use an AR model to 
represent the speech spectrum Assuming that speeoh signals can be accurately 
modeled as the output of an all pole filter driven by white noise for unvoiced 
speech, or driven by tram of impulses for voiced speech, the speech waveform may 
be reduced to a small set of parameters Thus, only the model parameters and the 
period of the impulse tram need be transmitted or stored Speech synthesis is 
achieved using appropriate model for each speech sound 

Note that we have considered the case where the sampling grid is not fixed 
but is arbitrary In a fixed sampling grid strategy, a globally optimum sampling 
grid is searched for and is used to sample the given signal In contrast to this, Me 
have considered the case where the sampling grid structure changes with the 
signal, depending upon the efficiency and/or noise reduction considerations 

Distributed Radar Deteotion Systems (DRDS) C73 make use of a fixed 
nonuniform sampling grid strategy The idea is to find globally optimum struotures 
and strategies for various elements and then logically combine all the respective 
signals in the so called central processor It is observed that the optimum 
peripheral detector thresholds depend upon the central processor strategies, 
(AND, OR, etc ) and vice versa 



CHAPTER 2_ 

FQRMULATPN OF THE PROBLEM AND MEASUREMENT 

INTERPOLATION 


z± jmBQDuenm - 


In this chapter we discuss the problem of transforming a nonumformly 
sampled data to uniformly sampled one, using orthogonal polynomial approximation 
This step will be necessary in order to be able to apply the parameter estimation 
techniques given m the subsequent chapters We also discuss, the effect of noise 
embedded in the nonumformly sampled data, on the approximating signal 

We have considered the more general case of complex exponentials It may, 
however, be noted that the measurement interpolation procedure itself is 
independent of the actual signal used Lastly we will discuss the relationship 
between the noise m the approximating signal to that of the nonumformly sampled 
signal 

22 FORMULATION - 


A complex signal r(t) which consists of M complex exponentials, can 
be expressed as, 


r<t) 



i=l 



, tio 


<2 2 1) 


For a real signal, the expression becomes a sum of damped sinusoids and we 
will have one complex conjugate pair of roots m the s domain for each of the 
damped sinusoids 

Let y<t> be the signal which is equal to r(t) added with corrupting noise w(t> 
Thus , y(t) can be written as, 


y<t) = r(t> + w<t) 



+ wit) 


, tiO 


(2 2 2 ) 





Let Ct k k = 1, 2, K), be the time instants where the signal y(t) is sampled Note 

that these time instants need not be at uniform spacing Uniform sampling becomes 
a special case of the given problem Thus, the corrupted signal y(t), at the instants 
t k , will be given by, 


y(t k ) 
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= z 

1=1 
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(2 2 3) 


We have assumed wttj,,) to be a complex, zero mean, Gaussian stationary 
sequence with known normalised autocovamance function White noise assumption 
for w(t k ) will be a special case where the normalised covariance function will be 
unity for lag zero, but will be zero elsewhere Also, as w(t k ) is assumed to be a 
stationary process, its normalised autocorrelation matrix will be symmetric and 
have a Toeplitz structure 

Thus, our problem is to estimate the parameters of r(t) using the corrupted 
samples y(t k ) 


23 ORTHOGONAL POLYNOMIAL APPROXIMATION - 


Our first 30 b will be to transform the nonuniform samples y(t k ), to a signal 
x(t), which will approximate y(t) over the sampling span If the sampling theorem is 
satisfied locally Eli than complete reconstruction of the original signal is 
possible We will be using orthogonal polynomials as given m Ell to generate the 
approximating signal x(t) Thus, x(t) will be given by, 

N 

x(t) = d 3 p jtt), <2 3 D 

3 = i 

where N is the order of the approximating polynomial 

For any arbitrary set of sampling grid Ct k k = 1, 2, K), we can evaluate 

the orthogonal polynomials p 5 from the recurrence relations given by, 

Pj - (t - V Pj-l “ 
with p 0 = 1, p. A = 0, 


b j-i Pj-2' 


(2 3 2) 
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y t k [ Pj.ivj 2 , 
k = i 
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k = i 


It should be noted that the this orthogonal set of polynomials mil be 
independent of the actual value of y(t k ) but will be dependent only on the 
nonuniform sampling grid Once the orthogonal polynomials p ? are calculated, we 
can also find d } To calculate 6 } we make use of equation <2 3 1) by replacing x(t> in 
it, by y(t), at t k 

Therefore, 

N 

*<V = 2 0,., V <2 3 3) 

3 = 1 


Rewriting equation (2 3 3) in vector form, we have, 
Y = PD 


(2 3 4) 


where, Y is the input vector given by, 
Y « [ y(t i ) y<t 2 ) 


P is given by. 


PqK^) p^tj) 

p 0 (t 2 ) P i (t 2 > 


P = 


P 0 (t K ) p^) 


y(tj<> 3 T 


p N-i (t l ) 

p N-i (t 2 ) 


(2 3 5) 


(2 3 6) 


p N-i (t K ) 



and D is the coefficient vector, given byj 


D - [ d L d 2 d N 3 T (2 3 7) 

Equation <2 3 4) can be solved for D, as both Y and P are known 

We now elucidate how the orthogonality property of the approximating 
polynomials can be used to find D Multiplying both sides of equation <2 3 4) by P T , 
we get, 

P T Y = P T PD (2 3 8) 

Equation (2 3 8) can be expressed as, 
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(2 3 9) 


T 

From equation (2 3 9), we observe that P P is a diagonal matrix This is 
because of the orthogonality property of the approximating polynomials The 
orthogonality property states that, 


2 p iV p a ( Hc > = °' 


for i # o 


= nonzero, for l = j 


(2 3 10> 


T -1 

We now multiply both sides of equation (2 3 8) by (P P) 


Thus, we get, 



D - <P T P) _i P T Y 


<2 3 ii) 


T T 

As P P as a diagonal matrix, it is evident that the inverse of P P is also a 

diagonal matrix, where each of the diagonal element is reciprocal of the 

T 

corresponding element of P P Thus, the choice of orthogonal polynomials effeots 
in reducing the computational burden required to solve equation (2 3 ii) 

Once dj are found, these can be substituted m the equation (2 3 i) to find x<t) 


Let c<t t< > be the difference between y(t k > and x<t k ) 
Therefore, 

c<t k ) = tKtjj) - X(t( < ) 


(2 3 12) 


Let e(t k ) be the approximating error Thus, e<t k ) will be the difference 
between x(t k ) and r(t k ) 

Thus, we get, 

e(t k ) = x(t k ) - r(t k ) <23 13) 


Therefore, 


x(ti,) = r<ti,) + e<ti,) 


<2 3 14) 


We observe, from equation <2 3 14), that the approximating signal x<t k ) 
consists of the original signal r<t k ) and the approximating error e(t k ) Also, from 
equation <2 3 12), it follows that, 

y(t k ) a x<t k > + ottfc) <2 3 IS) 

Substituting x(t k ) from equation (2 3 14) in equation (2 3 15), we get, 

y<t k > = r<t k ) + l c<t k ) + e<t k ) 3 <2 3 16) 

Comparing equations (2 2 3) and <2 3 16), we get, 

w(t k ) = t c<t k ) 3 + C e(t k ) 1 <2 3 17) 

Note that the SNR, as mentioned here, is the peak SNR This means that the 
SNR is calculated by considering the peak of the signal r(t) It follows that the 
local SNR will go on decreasing as the signal damps out 


* 


f 



S AMPLE 


r(t) is given by, 

S,t 5 2 t Sgt S 4 t 

r<t) = a^e + a 2 e + a 3 e + a 4 e , 
where, a^ = a 2 = a 3 = a 4 = 1 5, 

s L - -0 01322 + } 2*0 04, s 2 = -0 01322 - j 2*0 04, 

s 3 = -0 017664 + j 2*0 05456, s 4 = -0 017664 - j 2*0 05456 

SNR = 5 dB 


Fig 2 1 shows the plot for r(t k ) and x(t k ) Note that the samples for r(t k ) are 
nonuniform and arbitrary to the extent that the maximum inter sample distance 
does not fall below the minimum local Nyciuist rate required 

Fig 2 2 shows the plot for w<t k ) and e(t k > A very important oonolusion can 
be drawn from this graph Looking at the nature of e(t k ) we can see that it is not 
white, but is bandlimited 

From Fig 2 3 we see that c(t k ) also is not white, but is bandlimited Thus, by 
orthogonal polynomial approximation, we have effectively split the spectrum of 
w(t k ) m a upper frequency portion, which is c(t k ), and a lower frequency portion, 
which is e(t k ) The frequency at which this separation takes place is dependent on 
the polynomal order selected It is interesting to see from Fig 2 2 that e(t k > is 
approximately equal to the average of w(t k >, where the averaging process takes 
into account only a limited number of samples around a particular point 
Apparently, this number depends on the polynomial order seleoted Fig 2 2 is for 
polynomial order (N) of 22 

Fig 2 4 shows the graph of w(t k ) and e(t k ), but with N = 19 It oan be easily 
seen from Fig 2 4 that e<t k > is smooth as compared to e(t k ) in Fig 2 2 
Fig 2 3 and Fig 2 5 show c(t k ) for N = 22 and N = 19 respectively 

From Fig 2 2 and Fig 2 4, we can conclude that e(t k ) is present because of 
the low frequency components m w(t k ), which in turn introduces error in x(t k ) To 
support this statement we copy all the data points of c(t k ) in w<t k ) so that w(t k ) 
becomes bandlimited Now we perform the polynomial approximation procedure with 
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SNR of 5 dB Fig 2 6 shows the plot of r(t k > and x(t k ) in this case Note that r<t k ) 
has remained same but x<t k > has changed drastically as compared to the initial case 
of white noise, and is almost equal to Kt^) This is evident m Fig 2? which shows 
the plot of w(t k ) and e(t k ) From Fig 2 7 we see that e(t k > has a small variation 
compared to the variation m Fig 2 4 This is because of the absence of the low 
frequency components m w(t k ) 

Fig 2 8 shows the plot for r(t k ) and x(t (({ ) for the same case but with SNR of 
40 dB From these plots, we can see that the approximating error is negligible 


2 4 SELECTION CRITERIA FOR THE ORDER OF THE 



We will use the error variance as a measure for selecting the model 
order The error varianoe is defined C13 as. 



(2 4 1) 


where K is the number of nonuniformly sampled points and N is the order of the 
approximating polynomial 


There can be two criteria for selecting the order of the approximating 
polynomial 

(i) Select that order which gives a minimum error variance 
(n) Select that order after which the error variance does 
not decrease appreciably 

It was observed from simulation that if r(tj e ) itself has rapid variations in 
it, then both of the above mentioned criteria select the same order But if r(t k ) is 
smooth enough then the order selected through the first criterion is generally 
large as compared to the order selected through the second criterion It is also 
observed, that the amount of increase m error variance, is not appreciable, 
although the model order is reduced, while using the second criterion Thus, to 
save further computational burden, one should opt for the second criterion 






Fig 211 

SfcB = 40 dB. (Mhite noise.) 


It should be noted that, to find the correct model order, we have to find the 
error variance for all the model orders, starting from 1, to the maximum value 
possible (which is decided by the machine capacity) After this, select the 
appropriate order, and repeat the approximation procedure for that order only 

Fig 2 9 shows the plot of error variance vs polynomial order for a 
sequance w(t k ), whioh is white and with SNR of 5 dB Fig 2 10 shows the plot of 
error variance vs polynomial order for a sequence w(t^), which is colored and with 
SNR of 5 dB Fig 2 11 shows the plot of error variance vs polynomial order for a 
sequence wtt^), which is colored and with SNR of 40 dB 

25 OBSERVATIONS FROM SIMULATION - 


<i> For a given SNR, the nature of error variance vs polynomial order plot 
remained same for given rCt^) and for different w(t k ) sequences of 
same statistical properties, as expected The nature and looation of 
the oscillations in the error variance vs polynomial order plot 
remained approximately same 

(ii) The order of approximation, the nature of the nonuniform sampling grid 
and the nature of w(t k ), are the only factors that affect the nature of 
e(t k ) 

Note that the final order of approximation is decided by the nature of r(t k ) 
as well Thus, we can say that r(t k ) affects the nature of e(t k ) indirectly 

2£ REASON FOR THE PRESENCE OF THE APPROXIMATING ERROR - 


From the above discussion, it is clear that the orthogonal polynomial 
approximation procedure, seperates the high-frequency components, from the 
low-frequency components, present in the nonumformly sampled signal The 
frequency at which this separation occurs is decided by the polynomial order 
used For a more general case, where the nonumformly sampled signal is 
corrupted with colored noise, if the frequency spectra of the colored nosie and 
the actual signal, overlap m such a way that the colored noise spectrum goes 
^beyond the signal spectrum, then all those low-frequency components contained in 



the colored noise sequence, will appear m the approximating signal In case, these 
spectra are nonoverlapmg, then the error m the approximating signal is 
negligible 


27 THE IDEA OF PREPROCESSING - 


A very important observation should be mentioned at this stage The method 
provided m this thesis was found to give accurate estimates even for SNR levels 
as low as 7 dB This high accuracy can be attributed to the fact that the 
nonumformly sampled data was preprocessed (filtered) This, effectively reduced 
the noise present in the data, which was subsequently used by the parameter 
estimator Note that these preprocessed data samples have high SNR as compared 
to the SNR of the nonumformly sampled signal As a special case, the technique 
presented m this thesis oan also be applied for uniformly sampled signals having 
very low SNR, and very accurate parameter estimates can be found 

z&mmm - 


From the discussion in this chapter it is evident that the presence of 
error m the approximating signal, viz , e(t>, explains the inaccuracy in estimating 
the poles as shown by Sirkar and Sarkar 111, for considerable amount of noise 
Thus, if we can characterise this colored noise, we oan use this additional 
information to obtain more accurate estimate of poles, even when the noise level is 
quite high 





CHAPTER 3_ 

CHARACTER1ZATIQN OF ThE APPRQXIvlATING 

ERROR 


ai int roduc t i on - 


In the previous chapter we explained how a closed-form expression for a 
signal can be found by orthogonal polynomial approximation We also looked at the 
reason for the presence of the approximating error, and tried to illustrate that 
by simulation 

In this chapter we will move one step further by trying to find out the 
characterization of the approximating error, and the use of this information in 
parameter estimation 

It should be noted that all the simulation was done with white noise 
assumption for w(t k ) We will show in this chapter that the white noise assumption 
is not stringent and the derivations given in the subsequent sections are indeed 
generalizations for the oolored noise sequence also 

We further elucidate that if the so called normalized covariance function of 
the noise, embedded in the nonuniformly sampled signal, is known a prion, then the 
normalized covariance function of the error, m the approximating signal, can be 
estimated 

32 MATHEMATICAL CHARACTERIZATION OF THE APPROXIMATING ERROR - 


We rewrite equations (2 3 i) and (2 3 ii), respectively 
x<t) = 22 d J p J-i a) 

3 « 1 

D = (pWy 


(2 3 1) 


<2 3 ii) 


t 



Note that, p"*"p is a diagonal matrix, as is evident from equation (2 3 9) This 
is because of the fact that the approximating polynomials, p 3 , are orthogonal to 
each other (equation 2 3 10) Furthermore, (P^P)’ 1 will also be a diagonal matrix 
Thus, the computation of <P^P)’^ is simplified beoause of this orthogonality 
property 

Following conclusions can be drawn from equation (2 3 il), for a given 
nonuniform sampling grid and polynomial order 

<i) From previous chapter we know that P is dependent only on the 

nonuniform sampling grid and the polynomial order used 
T -i T 

<n) Thus, (P P> P will be a constant N by K matrix 
(in) Also, the coefficient vector D is solely decided by the input vector Y 

Substituting the vector D from equation (2 3 ii) m equation <2 3 1) we can 
s e g that we are effectively transforming the stochastio sequence y(t) to another 
stochastic sequence x(t) 

Let nT be the instants where x(t) is sampled uniformly (T is the inter sample 
distance) Therefore, equation (2 3 1) can be written as, 

xCnTl - ^ cljPj-ilnTl (3 21) 

j = 1 

For a given T, we can omit the notation of nT and replace it by n only We 
will adhere to this notation unless otherwise specified Thus, equation (3 2 1) oan 
be expressed as, 

A 

xCn3 = \ dj Pj.jM (3 2 2) 

3 - i 

Rewriting equation (3 2 2) in vector form, we have, 

X = P'D (3 2 3) 

where X is given by, 

X = C xCOl xCll xt23 xCL-13 ] T (3 2 4) 


and L is the total number of uniformly sampled points of x(t) 



hematical simplicity, we assume that. 


t 0 - 0 

to z nT < t^, 

< is the total number of nonuniformly sampled points of y(t) 


L by N matrix given by, 


PqEDI 

PiCOl 

P 2 C03 

PN-i C0] 

P 0 E13 

Pl [i3 

p 2 C13 

P N _iti3 

P 0 CL-13 

p^L-13 

P 2 EL-13 



<3 2 5) 


> given by equation <2 3 7) 

ltuting the value of D from equation <2 3 ii) m equation <3 2 3), we get, 

X * P' <P T P) _i P T Y <3 2 6) 


write equation (2 3 5) for Y 

Y = C ytt^) y<t 2 > y<t K ) ] T <2 3 5) 


It is evident from equation <2 3 5) and equation (3 2 6) that, we essentially 
form the nonuniformly sampled input vector, Y, to a uniformly sampled output 
r, X 


T -IT 

We know that (P P) P is constant for a given nonuni form sampling grid 

tit 

aolynomial order Thus, P' <P P) P 1 will also be constant Using this fact, we 
*>ee from equation (3 2 6) that the random vector Y is transformed to another 
om vector X by multiplying Y with a constant L by K matrix T, given by, 


T ■ P' (P T P) -1 P T 


<3 2 7) 


Thus, each of the random variables x[Q3, xCil, 
^nation of Y, given by the transformation, 


; k - 1 

= Tj qnm y<Wl > ' m = k - i 
^ m = 0 


, xEL-13 will be a linear 


<3 2 8) 


n, k-1 


the new variable, m, is introduced to avoid the notation of q. 


Note that n corresponds to the (n + i)th random variable at the instant nT 


The constant transforming coefficients q n m are such that, q i0 is the element 
» ith row and jth column of the L by K transformation matrix T 

From equation <2 3 9), (3 2 5) and (3 2 6), we can write the expression for q n mi 


Appendix A) 


Snm 



K <2 

£ ty 

3 = 1 


(3 2 9) 


where p jL . i En] s P^EnTl 


We know that the polynomials p 3 are deterministic Thus, Pj(t n ) will be 
elated with p 3 (t n+k ) 


PROBABILITY DISTRIBUTION FUNCTION OF THE APPROXIMATING 
ERROR - 


A probability distribution function (PDF) that is frequently used to model the 

istical behaviour of a random variable is the Gaussian Distribution We assume 

to be complex, white or colored, Gaussian, stationary stochastic process Thus, 

trtiplex random variable g(t k ) at the t^th instant, with mean and 

k 


anoe crg t , is distributed according to a Gaussian or Normal distribution if 


PDF is given by, 


"V 





random variable, 

+ 


g t , is given 

K 


by, 


(3 3 1) 



i Ut and Vi are independent real random variables, distributed as, 

K K 


I ag tL 

% ~ N [“ u v _ r 


N 




u \’ T 


>o, 


Mg t k " E ( 3 tJ * E KJ + 3E KJ 


= lhj. + 3 lL\j. 

l k l k 


id, 


var<g lk > = E jl (u lfe + 3 v^) - Ui u ^ + 3 


= E lu t ■ ii u | 


+ £ 


1 V t„ " J 


k ) 




ie 3 omt PDF of the vector G, where 
G = t g(t i > g<t 2 > 


g(t^) 3, is given by, 


p(G) = 


" ( G - Mq Cgg ^ ( G - IIq ) 


(n) K det( Cgg ) 


(3 3 2) 


(here ( 1 q is the mean defined as, 
*i Q = E £ G } 


(3 3 3) 


lS H 

(nd ( ) n is the operator for the conjugate transpose of a matrix 


the expectation operator, and C gg is the covariance matrix of the stochastic 
^rence G, defined as, 

C gg = E ( C G - tl G 3 C G - Uq 3 H ) 


(3 3 4) 



Note that Mg is a scalar whereas Mq is a vector 




We assume that the random sequence has zero mean Therefore, 
Mq = 0, null vector, 


(3 3 5) 


and 

Cgg = E ( G G H ) 

Substituting equation <3 3 5) and (3 3 6) in equation (3 3 2), we get, 

f- ( G H Cgg " 1 G ) 

p(G) s — „ L - e 

(71/ detj^ Cgg j 


(3 3 6) 


(3 3 7) 


Note that C gg will be a symmetric Toephz matrix, since, we have assumed 
w(t k ) to be a stationary process 

We now consider the effect of the polynomial approximation procedure on 
w(t k ) Consider the case when r(t k ) is a zero valued sequence Thus, y(t k ) will be 
equal to w(t k ), as can be seen from (2 2 3) 

Therefore, 

y(t k ) = 0 + w(t k > (3 3 8) 

Also, from equation (2 3 14) we can write, 

xCnl = 0 + etnl (3 3 9) 

Substituting equation (3 3 8) into equation (3 2 8) we get, 

K - i 

xCnl = ^ q n m w(t m+i ) (3 3 10) 

rn = 0 

From equation (3 3 9) and equation (3 3 10), we get, 

K-i 

eCnl = 23 q n m ^m+i) 
m = 0 


<3 3 11) 


From equation (3 3 11) we can see that the random variables etnD are not 
independent of each other This is because of the fact that the polynomials p 3 are 
deterministic Thus, Pj(t n ) will be correlated with P/t n+k ) Note that the 
orthogonality property does not come into picture, because we are considering only 
one polynomial, p 0 , at a time Thus, even in the case when w(t k )'s are uncorrelated, 
the approximating error sequence eEnl will be correlated This way, the white noise 
assumption for w(t k ), then, becomes a special case 


In equation (3 3 ii), we have used the same polynomial order K, as was used 
when r(t k ) is a nonzero sequence Thus, if w(t k > are known, then eCnl, which is the 
approximating error, can be found accurately But this statement is eluding in the 
sense that, if w<t k ) was known then we would have subtracted it from y(t k ) right at 
the start Thus, it becomes necessary to find the statistical properties of eCnD for 
parameter estimation We will now try to investigate this possibility 


Equation (3 3 7) gives the joint PDF of the random sequence w<t k ) We can find 
a similar expression for PDF of eCnl, as given by <3 3 ii) 

Note that the joint PDF of eCnl will also be Gaussian, because, from equation 
(3 3 ii) we can see that eCnl is a linear combination of the random variables w<t k ), 
and each of the random variables w(t k > is assumed to have the same 
Gaussian PDF [43 

Thus, the PDF of eCnl can be written as, 

f- ( E H Gee ' 1 E ) | 

p(E) = —= i- e J (3 3 i2) 

ta) L det( C ee ) 

where, 

E = [ e£03 eti3 eCL-i3 ] T <3 3 i3) 


and 


C eB = E ( E E H ) 


<3 3 15) 


3L1 EflE SE THE EBF Q R EBB NON-G4 /SSi,4N m.tt 


We will make use of the Central Limit Theorem to find the PDF of the 
approximating error, when w(t k ) has a non-Qaussian PDF 

The Central Limit Theorem states that, for the random variables uj, under 
general conditions, the density p(u) of their sum, 
u = + u 2 + + u n 

properly normalized, tends to a normal curve as n -t » 

In other words, if n is sufficiently large, then 

- C - u -^' 2 ) 

p(u) = — e Cru 

it ffu 


It has been shown E43 that the theorem does not hold if only a small number, 
m, of the given densities are dominant, le, if the densities of the other random 
variables are relatively narrow m the sense that, in the evaluation of the 
convolution, they can be approximated by impulses Furthermore, it has been stated 
[43 that p(u) is effectively the convolution of only the m dominant densities and it 
need not be close to a normal curve However, if the densities p(u JL ) are smooth, 
then p(u) is nearly Gaussian even when n is small 

In our problem, the value which n takes is decided by the polynomial order N 
used for approximation Furthermore, this order depends upon the amount of 
variations in the input signal 

From the simulation results given m the subsequent chapters it is shown 
that the Gaussian assumption for eEnl indeed holds true, for all practical 
purposes 



3J5 THE NORMALIZED CW/iRI/lNCE / CORRELATION MATRIX 


Following are some of the properties of the covariance matirix 

(l) If the stochastic sequence is stationary and white, then the covariance 
matrix is a diagonal matrix with equal elements The diagonal elements 
are equal to the variance of the process Note that, as the process is 
assumed to be stationary, each of the random variables will have equal 
variance 

(n) If the stochastic sequence is stationary and colored, then it's covariance 
matrix is a symmetric Toeplitz matrix 

(in) For a white nonstationary stochastic sequence, the covariance matrix 
is a diagonal matrix but with unequal elements 
(iv) For a colored nonstationary stochastic sequence, the covariance matrix 
does not posses any of the special properties mentioned above 


We define a nonstationary process as the one which is not a Wide Sense 
Stationary (WSS) Ml process A process is called WSS if its mean is constant and 
its covariance depends only on the lag parameter and not on the exact time 
instant 

We now define the so called Normalized Covariance Matrix Note that the 
Normalized Covariance Matrix can be defined only for a WSS process 

The Normalized Covanane Matrix is that matrix which is equal the matrix 
formed by normalizing all the elements of the Covariance Matrix with respect to 
any of its principal diagonal elements , all the principal diagonal elements of the 
Covariance Matrix will be equal, for a stationary process In other words, we get 
the Normalized Covariance Matrix by dividing each of the elements of the 
Covariance Matrix by the variance (for a zero mean process) Note that the 
covariance manx reduces to a correlation matrix for a zero mean process 

The normalized covariance matrix is an identity matrix for a white noise 
stationary sequence It will be a symmetric Toeplitz matrix , for a colored noise 
stationary sequence 

We can express equations (3 3 7) and (3 3 12), respectively, m terms of their 
normalized covariance matrices and variances 







< 


Therefore, 


and, 


p«3> = 


p<E) = 


OjCtFg ) Cfst.| J 


Cittre 2 ) 1 * det| C' ee } 


- ' Q H C'gg " 1 Q ) 

■ 




Aj < E H GW 1 E ) 


<3 5 i) 


<3 5 2) 


where, C'gg and C J ee are the normalized covariance matrices for g<t k ) and etnl 
respectively 

In Section 3 3 we stressed on the need to find the statistical properties of 
eCnl We will now elucidate how this can be done 

Because of zero mean assumption, the covariance function is equal to the 
correlation function The correlation matrix R gg for w(t) can be written as, 


rggCOl PggjH.1 


rggC-Kl 


R gg 


rggtil PgglOl 


r g g[-(K-i)l 


(3 5 3) 


rggCKl rggCK-il 


rggCOl 


where K is the number of uniformly sampled points 


But as g<t> is WSS, equation <3 5 3) can be written as 


rggC03 rggCil rggCK] 

r gg Ci3 PggCOl r g gCK-il 


rggflO PggCK-i3 


r gg C03 


since, 

r gg*tk3 — pggC-kl , 

where, f denotes the conjugation operator 


(3 5 4) 


(3 5 5) 



Similarly, the correlation matrix R aa can be written as, 


reglOl r ee [13 


r aa [L3 


r aa C13 r ea C03 


'ee = 


r ae [L3 r ee CL- 13 


r ee i:L-13 


r ae [03 


(3 5 6) 


Note that, both R gg and R ea are symmetric Toeplitz matrices, and r gg Ck3 and 
r e otk3 are the respective correlation functions, at lag k 

We now rewrite equations (3 3 7), (3 3 12), (3 5 1), (3 5 2), by substituting the 
correlation matrix R, instead of the covariance matrix C 
Therefore, we get, 


p(G) 


- < G H Rgg" 1 G ) 


(•n) K det[ Rgg ] 


p(E) = 


( E H R ee _i E ) 


(tt) L det( R ee J 


p(G) = 


- — ( G R R’gg ^ G ) 1 

J 


(TTCFg 2 r det( R'gg ) 


p(E) = 


1 

(it<7 B 2 ) L det[ R' ee J 


e 


-K ( eH R ee _i ^ > 
er e 


(3 5 7) 


(3 5 8) 


(3 5 9) 


(3 5 10) 


Note that R’gg and R’ e e are the normalized correlation matrices of gd^) and 

eCnl respectively Also, 

2 2 
og = r gg I03 and cr e 


r ee C03 


Following relations hold for the WSS processes g[n3 and eCn3 respectively 


The mean is given by, 


/lg = £ 1 

( gCnl j , 

(3 5 11) 

{Iq = E I 

[ etnl ) 

<3 5 12) 

The autocorrelation functions (ACF) are given by, 


rggtkl = 

E ( g*Cn+k3 gtn3 j 

(3 5 13) 

r ee tkl = 

E ( e*tn+k3 etn3 ) 

(3 5 14) 

The crosscorrelation function (CCF) is given by. 


rgetkl = 

E ( g*Cn+k3 etn3 ) 

(3 5 15) 


For a linear shift invariant (LSI) system with impulse response hCnl and with 
a WSS random input, various relationships between the input and output correlation 
functions and htnl can be given 


Let g[n3 be the input sequence and eCnl be the output sequence 
Therefore, 


r ge Ek3 = htkl ® r gg [k3 


oo 


= hCl] rggCk-11 

l T%o 


regChl = h*C-k) ® rggtkl 


(3 5 16) 


r e eCW = htkl © r a gCk3 


= hCk3 @ h*C-k3 ® r gg £k3 
oo oo 

= 22 h[k-m3 22 h*[-13 r gg Cm-13 (3 5 IS) 

ms^oo i 

where © denotes the discrete convolution operator 

We now consider how to find the correlation function, r ee Ck3, of the 
approximating error e[n3 

* 

First step is to find h C-k3 and h£k3 from equation <3.2 6) Next we find r ew Ek3 
from equation <3 5 17 ) Now we can use r eg Ck] to find r ee Ck3 from equation (3 5 18) 
Once r ae Ek3 is known, it can be used to find the correlation matrix R ee from 
equation (3 5 6) 

Note that the computaion of R e e as Siven above will necessitate us to know 
the exact value of r ww Ck3 in the first place This may not be possible m many 
cases This problem is solved by using the normalized correlation function instead 
of the actual correlation function of wCn3 We denote this normalized correlation 
funotion of wtn3 by r' ww tk3, defined by, 

<3519> 

Similarly we define the normalized correlation function, r eaEkl. of eCnl by, 

r' ee Ck3 = (3 5 20) 

ee r ee t03 

Note that this so called normalized correlation funotion preserves all the 
statistical properties of any WSS process 

The normalized correlation function can also be interpreted as the 
correlation function of the same stoohastio process which has unit variance 

Furthermore, it is evident that if wCn3 is white then its normalized 
correlation function will be such that, 


i 




r'C03 

r’Ek3 


i, and 

0, for all k & 0 


1 


and the corresponding normalised correlation matrix R' will be an identity matrix 

Similar interpretation oan also be made for any colored WSS process 

From the above discussion we observe that for a WSS process, rCOl is 
proportional to the energy content in the process Furthermore, any change m the 
energy content will change the variance of the process leaving the normalised 
correlation function unaffected 

In subsequent chapters we will develop methods which will use the 
normalized correlation function estimate, as given by equation (3 5 20), to get 
accurate parameter estimates 



CHAPTER 4. 

PARAMETER E8TIMATDN IN COLORED GAUSSIAN NOISE 


4 A INTRODUCTION - 


In the previous chapters we studied the transformation from the nonuniform 
sampling and domain, to the uniform sampling grid domain We also saw how colored 
noise is introduced because of the polynomial approximation procedure We also 
proved that, if the signal and noise frequency spectra do not overlap, then this 
colored noise is negligible The problem of finding the statistical properties of 
the colored noise m the approximating signal was discussed fully m the previous 
chapter 

The purpose of this chapter will be to use the statistical properties of the 
colored noise during parameter estimation The signal is modelled by a sum of 
complex exponentials, embedded m colored noise whose normalized covariance 
matrix is known All the simulation was done on real data, which essentially 
becomes a special case 

Lot of pioneering work in parameter estimation for complex exponentials has 
been done by Roy and Kailath C61, Kumaresan and Tufts C93, and Farrier et al CiOl 
The technique discussed here considers the more general case of colored noise 
Furthermore, this technique exploits the properties of the noise corrupted signal 
autocorrelation-like matrix by using subspace separation approach 

It is well known that a signal consisting of complex exponentials can be 
modelled by a difference equation We have presented a technique which uses this 
difference equation, and achieves accurate estimates of complex exponentials using 
Singular Value Decomposition till 

Further enhancement m accuracy is achieved by using the signal subspace 
seperation approach, which is also referred as Principal Components Analysis 1121 


±2 tmj^M UC AL M WJCTf-mtriOtt QE WK JPPROXiMmM signal 


The approximating signal, x<t), consisting of M complex exponentials, can be 
expressed as. 


x(t) 



e<t) 


(4 2 i) 


Let xCOl, xtil, , xCL-i] be the uniformly sampled corrupted sequenoe 

Thus, the discrete-time prooess, xCn], can be modelled by a time-senes as, 


xCnl = rCnl + etnl 

where, rCnl is the nth sample of the nonoorrupted signal term, given by. 


rCnl 



alkl rtn - kl 


<4 2 2) 


which is the difference equation representation of the Prony's model Ell The 
coefficients aEkl are such that the impulse response of the filter, Ate), given by. 


Ate) - 


1 + aC13z"* + aC23z"^ + 


+ aLMlz 


-M 


(4 2 3) 


s.n 

is equal to > Aj e , for n i M 


k = i 

Let s 1 and Z* be the roots in the s and z domain, respectively The transformation 
from the continuous to discrete-time domain is achieved through, 


e 


s l T 


<4 2 4) 


where, T is the inter sample time period, and s^ is given by, 


Sj s + jWj = « i + j2nfj (4 2 5) 

where, a l is the damping factor and f 4 is the normalized frequency 

The error, eCnl, is an innate part of the model and gives rise to the random nature 
of the observed process xCnl In our case, this is the colored approximating error 
sequence 




From equation (4 2 2), it follows that, 


xCM3 

= -aC13 rCM-13 - 

- aCM3 rt03 

+ 

eCM3 

xCM+13 

= -aC13 r£M3 - 

- aCM3 rE13 

+ 

etM+13 

xCL-13 

= -aC13rtL-13 - 

- aCM3rCL-M-13 

+ 

eCM-13 


Equation (4 2 6) can be written m vector notation as. 


X = QA + E 


(4 2?) 


where, 


X = [ x[M3 xtM+13 


xEL-13 ]' 


(4 2 8) 


Q = 


rCM-13 

rCMD 


rCL-23 


K03 

rC13 


rEL-M-13 


(4 2 9) 


and 


A = [ -at 13 -at23 


E ■ [ etH3 eCM+13 


T 


-aCM3 ] 


eCL-13 ]' 


(4 2 10) 


<4 2 11) 


The joint PDF of E, as given by equation (4 2 11), can be written as, 


p(E) 


(itcr e 2 ) L det (R' ee ) 


e 


-ij 


(4 2 12) 


' where, the normalized correlation matrix, R'ee' is found by the method mentioned in 

% 

* the previous chapter Substituting equation (4 2 ?) m equation (4 2 12), we can write, 


p(X, A) = 


( X - QA ) H R , ee " i ( X - QA )J 


(tt c e 2 ) L det (R' ee ) 


<4 2 13) 


where, the notation p(X, A) explicitly shows the dependence of the joint probability 

density function on both X and A Note that X, Q and R' ee are known auantities, 
2 

while A and «? e are unknowns 

We find the maximum) likelihood estimate (MLE) of A by maximizing equation 
<4 2 13) with respect to A C123 This can be done equivalently, by minimizing, 

S = ( X - QA ) H R'ee" 1 ( X - QA ) (4 2 14) 


which is a quadratic function of A After minimizing equation (4 2 14) with respect 
to A, we get the estimate of A as, 


A = ( Q H R'ee" 1 Q ) * Q H R 


H o, -i Y 
ee * 


(4 2 15) 


and the value of the quadratic function, S, is. 


S MIN = xH ( R'ee" 1 - R'ee" 1 Q ( Q H R'ee" 1 0 ) _i Q H R'ee" 1 ) X (4 2 16) 

We rewrite equation (4 2 15) as, 

( Q H R' e e -1 Gl ) A = Q H R'ee" 1 X (4 2 17) 


Therefore, 

T A = B (4 2 18) 

where, 

T = ( Q H R'ee" 1 Q ) (4 2 19) 

and 

B = Q H R'ee" 1 X (4 2 20) 


Thus, r can be interpreted as the estimate of the signal autooorrelation- 
like matrix, where the signal is corrupted by colored noise 

We can write the noise corrupted signal autocorrelation-hke matrix, F, as, 
(See Appendix B) 


r 


Rss + Ree 


<4 2 21) 



We perform spectral decomposition of R ss and R' BB Thus, we get, 


r = X x x 5 u i u i H + ** 2 x i e u i u i H <4 2 22) 

i=l i=l 

where, X 1 s is the ith eigen value of the signal autocorrelation-like matrix and X a e 
is the ith eigen value of the colored noise normalised autocorrelation matrix, and 
both R aa and Rss are of dimension M by M, 11 being the number of complex 
exponentials present m the corrupted signal Note that both the signal and noise 
subspaces are spanned by the same eigen vectors, u t , of T This is because an 
M-dimensional space can be spanned by any set of M orthogonal vector in the 
space We observe that, as the signal is locally WSS [43, the column and row eigen 
vectors of R ss are same (See Appendix C) Same argument is also applicable for 
the autocorrelation matrix of colored noise 


In almost all practical applications, the value of M is unknown Thus, we 
assume some large value, p, possibly much larger than M 
We rewrite equation <4 2 22) as, 


- X 

i = i 


u i u i 


H 


+ cr e 


P 

X 

0 = l 


U 1 U 1 


H 


(4 2 23) 


As R ss is of dimension p by p and it's rank is only M <H « p), it is not full 
rank The sample autocorrelation-like matrix, T, is, however, full rank due to the 

2 © g 

inclusion of the cr e term In the special case, when etnl is white, ^ are unity 

The frequeny and damping factor information about the signal is contained 
m the eigenvectors Since, r, as given by equation (4 2 19), is not the actual 

A 

signal autocorrelation matrix but an estimate of it, the estimate, r, can be 
expressed as, 



A 2 



A 

V, 



+ 




A 

V, 



(4 2 24) 


This is evident from the fact that, 

X A S = 0, for i = M+l, , p 

It should be noted that, r, as defined by equation (4 2 19), is not the exact 
value of the noise corrupted signal autocorrelation matrix To avoid the 
discrepancy to propogate any further, we have used p\, s (instead of X, s > m the 



principal eigen values in equation (4 2 24) Furthermore, the eigen vectors, v 1 , are 
related to the signal eigen vectors, e A , by the equation, 


A 

V 



p 

“1 


(4 2 25) 


4 3 SINGULAR VALUE DECOMPOSITION - 


We will use SVD for three reasons 

(i) To form the spectral decomposition as given by equation (4 2 24), 

A 

and use it to find the estimate, A 
(n) To find the model order M, from the spectral decomposition 
(n) To seperate the signal and noise subspaces and thereby achieve 
further accuracy m parameter estimates 

By applying SVD Cl 13, we decompose a n by m matrix A, such that, 


A = U V H 



>.4 3 1) 

where, 




U = t u i u 2 

U 1 u l+l 

u n 3 

(4 3 2) 

V — [ Vjl v 2 

V 1 v l+l 

v n 3 

(4 3 3) 

and 




= diag< Sjl S 2 

$!_! $1 0 0 

0 ) 

(4 3 4) 

where, £ $ 2 £ 

£ ti > 0 , 




and 1 is the rank of the matrix A 

-1 H X LI 

Also, U and V are unitary matrices, i e , U" = U and V" = V , respectively 
The column and row eigen vectors, u JL and v i( respectively, are defined as, 



\ 1 


l_J 

'lultiplying equation <4 3 5) by A and using equation <4 3 6), we get, 

A H A u = £ A H v = £ 2 u (4 3 7) 

Thus, u is the eigen vector of A A, with $ as the corresponding eigen value 
"lultiplying equation (4 3 6) by A and using equation (4 3 5), we get, 

A A H v = 5 A u = £ 2 v (4 3 8) 


H 

Thus, v is the eigen vector of A A , with <? 


as the corresponding eigen value 


If A is a square matrix, then one can find its eigen value decomposition 
(EVD) We will see how SVD can be used to find the EVD 

LJ Lj LJ 

From the basic properties of EVD, we know that A (A ) and A A (AA ) have 

the same eigen vectors Also, if X is the eigen value corresponding to an eigen 

H H H 2 

vector of A (or A ), then the eigen value of A A (or AA ) is equal to X Thus, if A 

is a square matrix, then its eigen value is the square root of the corresponding 

singular value [ill 

A 

For the noise corrupted signal autocorrelation-like matrix, T (whioh is a 
square p by p matrix), we can use SVD to find the spectral decomposition as given 
by equation (4 2 24) Note that f is found using equation (4 2 19) Also, as f is 
symmetric, the row and column eigen vectors will be same 


Substituting the value of f from equation <4 2 24) in equation (4 2 17), we get, 
( X ( p V 3 + V ) °1 °1 H + 2 V °1 ° 1 H I ft ■ ° H R 'ee _i X 

^1 = 1 1 ss M+l * 


(4 3 9) 



Equation <4 3 3) can be written in matrix notation as, 



Therefore, 

[0 £ 0 H 3 A = Q H R'ee -1 X (4 3 ±i> 

where, 

0 - [ 0 X v 2 v M 0 P ] (4 3 12) 


and 



A A -1 A H 

As V is a unitary matrix, V = V 
Therefore, equation (4 3 11) becomes, 


A = ( 0 ~ L V H ) Q H R'ee' 1 X 


(4 3 13) 



(4 3 14) 




r~ j 


rherefore, 


A 


_ V v i v i 
1 = 1 1 


i A A H I nH d* "i 


Q R'i 


'jhere X x are the eigen values found from the EVD of f 


(4 3 15) 


Example - 


rCn] a 1 5 e 


- ( ^ + j 2 % 0 08 ) nT 


+ 15 


-(e|j + 02ttOii)nT 


+ 1 5 e 

K = 50, L = 100 


- ( ^ - 3 2 it 0 08 ) nT 


+ 1 5 e 


(sj-j2it0 1i) nT 


Fig 4 1 shows the pole estimates using equation 4 3 15, for 10 different 
trials The white noise corrupted signal is sampled at 50 nonumformly spaced time 
instants As the data consists of two sinusoids, we have two complex conjugate 
pairs of signal poles in the z and s domains The estimate of the normalised 
autocorrelation matrix, R' ee , is calculated using equations 3 5 i?, 3 5 18 and 2 3 19 
The SNR in this case is 40 dB and the value of the model order, p, is 10 Fig 4 2 
and Fig 4 3 show the pole estimates for SNR of 20 dB and ID dB, respectively, with 
pc 10 (Mention about the actual pole locations) 

Figures 4 4, 4 5 and 4 6, show the pole estimates for SNRs of 40 dB, 20 dB 
and 10 dB, with an increased value of p, equal to 16 

The table given below shows the performance of the estimator for different 
noise levels and with varying values of the seleoted model order, p The estimation 
procedure was repeated for 10 trials, thus, the bias and variance, mentioned below, 
will give a rough understanding about the statical performance of the estimator 


Table £1 - P ■ 10 


LZ 

SNR = 40 dB 

SNR = 20 dB 

SNR = 10 dB 

1 

Bias 

Variance 

Bias 

Variance 

Bias 

Variance 

6 234 x 10“ 4 

9 959 x 10" 7 

9 203 x 10" 4 

1 159 x 10" 7 

1 651 x 10" 3 

4 073 x iO" 5 

m 

8 03? x 10" 4 

5 126 x iO -8 

9 127 x 10" 4 

2 771 x 10’ 7 

1 84? x ±0" 3 
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7 574 x 10" 3 

1 574 x iO“ 3 
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la £2 - p = 16 


SNR 

= 40 dB 

SNR = 

20 dB 

SNR 

= 10 dB 

Bias 

Variance 

Bias 

Variance 

Bias 

Variance 

2 576 x IO" 4 

1 174 x 10" 7 

7 880 x 10" 4 

8 303 x 10" 6 

1 418 x 10' 3 

1 418 x iO" 5 

7 659 x 10" 6 

2 675 x 10 -9 

8 240 x IO -5 

2 360 x 10‘ 7 

ESBBm 

mmm 

4 118 x 10" S 

2 859 x 10” 7 

6 330 x 10" 5 

2 235 x 10‘ 6 

5 521 x 10' 4 

7 06? x iO" 6 

5 677 x 10" b 

5 772 x 10 -9 

1 492 x IO -4 

2 649 x IO" 7 

2 621 x 10' 4 

1 027 x 10 -6 


Following conclusions can be drawn from the location of the pole estimates, 
different SNR levels and different values of p 

i> The pole estimates are accurate down to the SNR level of 10 dB 
l) As the SNR is decreased, the variance of estimated poles 
increases accordingly These perturbations are introduced 
because of the ill-conditioning of the problem C93 In Section 4 5 
we will discuss how this ill-conditioning can be taken care of, with 
corresponding increase in accuracy 
i) The variance of the estimated poles is reduced by increasing the 
model order, p We will elucidate this fact m the subsequent section 

MODEL ORDER SELECTION THROUGH SPECTRA- DECOMPOSITION - 


In this section we will address the problem of model order selection from a 
>en data set Equation (4 3 10) will be used for this purpose From equation 
3 10), the principal eigen values, viz , C pX a s + o e 2 \ e i = 1, , M ), can 

selected, such that they are large as compared to the nonprincipal eigen 
lues, viz , C a e 2 X A e l = M+i, , p 5, through proper ohoioe of p Note 

it, as Xj 6 are the eigen values of the normalized autocorrelation matrix of the 
proximating error, they are very small in magnitude, when the noise variance, 
% is also small Thus, the M principal eigen values are distinctly large than the 
w M nonprincipal eigen values, for a proper choice of p This fact can be used 
^appropriately selecting the model order, M 


IF 





















ram equation (4 2 24), it is evident that the estimate of the autocorrelation-like 
latnx, f, oan be written as, 


AAA 

r - A + B <4 4 i) 

A A 

ihere, A and B are defined to be the principal and nonpnnoipal terms, 
espectively 

We maintain that principal and nonpnncial eigenvalues are easily 
listinguishable This fact is due to the fundamental property of spectral 
lecomposition, that, 

min || f - A || = X z+1 (4 4 2) 

A A 

or any z < p, where, X x are the eigen values of f 

A 

In another words, a full rank matrix T, whose rank is p, oan be approximated 

A 

>y a low rank marix A, whose rank is z, where z < p, and the approximation error is 
liven by equation (4 4 2) 

Thus M corresponds to that minimum order after which the rate of change of 
he Mth eigen value, with respect to z, is almost constant, i e , 

dCmmlir - AIM constanl 
dz 

Je define this rate of change as, 

d E min II r - A || ] _ ^z+1 ^z+2 ^ ^ gj 

dz “ i 


Thus, z is a correct order, if it corresponds to that minimum order, after 
tihich the difference between the (z+i)th and <z+2)th eigen value is small as 
compared to the difference between the zth and (z+i)th eigen values 

During the process of simulation, it was observed that equation <4 4 3) is not 

t 

/ary usable for model order selection One reason could be that the estimate of 

he autocorrelation-like matrix is not exact, due to limited amount of data With 

f 

increased amount of data, it is only expected that equation (4 4 3> may give 
j|»curate model order estimates 


An alternative approach for model order selection has been proposed by 
idzow and Wu [133 which examines the norm ratio, R The norm ratio, R, is defined 
}> 


R = 



+ Xg 2 
+ Xp 2 


i 


(4 4 4) 


nere, q is the actual model order, and p is the model order used for simulation, 
ith p > q, and CX A 1 = 1, 2, p3 are the corresponding eigen value estimates 

F the autocorrelation-like matrix We can easily see that the norm ratio R always 
es between zero and one Moreover, if M readily identifiable largest eigen value 
stimates exist, then the smallest value of q for which the rank ratio, R, tends to 
nity, will correspond to the actual model order, H Ci.33 


Practically it is observed that the the ratio of the largest principal 

mgular value, to the smallest principal singular value, of the estimated 

2 

utocorrelation-like matrix, is large (Note that X x = for i = 1, 2, , p ), 

or any SNR level The transition from to is evident only for moderate 

NR levels, of about 20 dB Selecting the appropriate model order becomes 
iffioult for SNR below this threshold, using traditional methods But as we have 
reprooessed the data by using orthogonal polynomial approximation, the SNR of 
he approximating signal is increased, comparatively Due to this fact (it was 
'erified from simulation), the transition from to £ M+i , xs svident even uptill 
>NR of 7dB, when the data is preprocessed Note that, m our case the polynomial 
ipproximation procedure was necessary, to enable the transformation from 
lonumform sampling grid to uniform sampling grid Furthermore, we have shown 
hrough simulation in Chapter 5 that, by using the idea of preprocessing, for a 
jniformly sampled data, we get very accurate pole estimates as compared to 
raditional methods 


Keeping the practical considerations m mind, we propose a slightly 
different criteria for selecting the appropriate model order 

We define a criterion C, such that, 


C = , for . = 1. 2, 


P-i 


(4 4 5) 


iUS, the order, i, will be equal to the actual model order, M, when C is maximum 
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Fig 4 7 shows a typical plot of the norm ratio, R, vs the model order, for 
SNR of 7dB, and with p equal to 20 <Note that the value of p is large) M was 4 for 
this simulation Fig 4 8 shows a typical plot of the criterion, C, vs the model order, 
for the same example Comparing these two cases, we observe that the actual model 
order, M, can be easily distinguished by using the criteria given by equation 
(4 4 5), for the example under consideration Fig 4 8 and Fig 10 show the plots for 
the R and C vs model order, respectively, but for p = 16 We see from Fig 4 9 and 
Fig 4 10 that it is more difficult to select the model order when p is reduced 

Note that the value of p should be seleoted in such a manner that the 
smallest principal eigen value (i e , the Mth eigen value) is sufficiently large than 
any of the nonprincipal eigen values This is to say that, any decrease m the SNR 
level should be accompanied by a corresponding increase in the value of p, so that 
the principal eigen values are easily distinguishable If this condition is not met, 
then model order selection becomes an impossible task, for very low SNR levels 

4.5 THE SUBSPACE SEPARATION APPROACH - 


Upto this point, we discussed how to find the parameter estimate, A, using 
equation <4 3 13) In the previous section we also saw how spectral decomposition 
helps us m selecting the model order In this section, we will formulate a 
technique, which achieves further accuracy in the parameter estimates, than that 
can be achieved using equation (4 3 13) 

We will assume that model order selection has already been done The order 
selected will correspond to the actual order, M, only if the SNR is not below, 
say 7 dB This statement has been explained in the previous section 


Now that M is known, we can effect the subspace separation, as given by 
equation (4 4 1), where A uses only the principal eigen vectors and eigen values As 
M corresponds to that minimum order for which the error norm, defined by 
equation <4 4 2), is minimum, we oan approximate f , by A < whioh is of rank M> Thus, 
the new estimator can be written, such that it uses only the M principal eigen 
vectors and eigen values Thus, 
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A similar approach (but restricted to the white noise case) has also been used in 
Ref 9 The subspace separation approach mentioned here is referred as truncated 
gj/jj m Ref 9 The principal difference between the approach given here and the 
approach used in Ref 9 is that, we use an autocorrelation-like matrix associated 
with the time series unlike the data matrix itself as done m Ref 9 It was 
observed that the autocorrelation-like matrix gives more accurate estimates as 
compared to the estimates from the data matrix itself This is due to the fact that 
the noise effects are considerably moderated by using the autocorrelation-like 
matrix along with a large value of p 

Fig 4 11 shows the pole estimates for 20 different realizations of the 
random sequence, w<t k >, using the estimator given by equation <4 5 1), with SNR of 7 
dB and P equal to 20 Fig 4 12 shows the pole estimates for 20 different 
realizations but using the estimator given by equation <4 3 15), for the same case 
Clearly, it oan be seen that the variance of the pole estimates shown m Fig 4 11 
has reduced substantially as compared to the variance of the pole estimates shown 
in Fig 4 12 We also observe a certain amount of radially Inward shift in the 
location of the pole estimates on the z domain of Fig 4 11 This phenomenon of 
radial bias will be explained in detail in section 4 7 

4 6 FURTHER TMPRIWEMENT S. IN UM SUBSPAQ E SEPARrfTIOJJ A PPROACH _ 


We write the spectral decomposition as given by equation (4 5 1) 


A 


r 



*e 2 ) 


A A 

v l v i 


H 


(4 61) 


where, the estimated eigen vectors, 0 t , are related to the actual signal eigen 
vectors, e^ by equation (4 2 25) See C143 for properties of the autocorrelation-like 

A 

matrix, r 


The matrix B, whioh is spanned by the p - M nonprmoipsl eigen veotors is 
given by, 

6 . V V s «i »i H <4sa 

i rrf+i 

Using the estimate R'ee' 1 * w * can effect lVs speclral decomposition, to find the 
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estimates ! i ~ L, 2., , p y Using these eigen values in 

equation <4 6 2), we can find the estimate of the colored noise variance, o e ^ It was 
observed during simulation that the noise variance estimate is correct, only when 
p is moderately high Thus, using the colored noise variance estimate, & e ^, and 
£ i = i, 2, , M ■», we can subract the terms C i = ±, 2, 

, M 3, from the M principal eigen values of T, as given by equation <4 6 1) 
Thus, for high values of p, it is only expected that the autooorrelation-like matrix, 
f , can be written as, 




r 




<4 6 3) 


Thus, the new relation for finding the parameter estimate, App^, is given by, 


A MPC 



/ / H 

v i v i 


Q H R'ee" 1 x 


<4 6 4) 


The pole estimates, using the estimator given by equation (4 6 4), for 20 different 
realizations, are shown m Fig 4 13 M is 4, p is 20 and SNR is 7 dB Comparing Fig 
4 11 and Fig 4 13, we observe that the t aclial bias has reduced The actual pole 
looations are shown by cross in Fig 4 9 and Fig 4 11 


47 RADIAL BIAS - 


We will explain the phenomenon of radial bias using equations (4 6 1) and 
(4 6 2) ( In Ref [93, the phenomenon of radial bias, is depicted, for white noise 
case, on page 837) 


U 

The eigen vectors, e A and e 1 n , of the signal autocorrelation matrix, Rgg, can 
be written as, 

T 


e i = 


and, 


2(a 1 +jo 1 ) 

1 e e 


(p-l)(« i +jaj 1 ) ] 


H 


Kfltj-jUj) 2(a i -ju 1 ) 

1 e e 


(p-lXa^-jt^) 


(4 7 1) 


and, 

r r 


M 


: ss = & L e x H 


i = l 


(4 7 2) 





where, 


\ - A i = Power of the ith sinusoid m the signal sEn3, given by, 


M 


and, 


sEn3 = 2l, a i e 

1 = 1 


s x *5 a x + 3 


(4 7 3) 


(4 7 4) 


From equation <4 6 i), we can see that the estimate of say X' 4 , which we 
calculate, is always greater than the actual value of \ 

Therefore, 

X i > \ 

Thus, for any ith sinusoid, 

» 2(e 


(4 7 5) 


X'i Ui Ui H = [ \ s + ^L. Lu” 

. j 


> \ U A U A H (4 7 6) 

From equation (4 7 6), we observe that, the increase m the estimated eigen values 

<7e 2 X e 

is due to the term — p Also, if eEn3 is a white noise sequence, then X^ e will be 

unity 

Thus, we can write, 


ie, 


s,n s.n 

\\ e > \ e 1 


a.n j 2-Kf.n a.n jZnf.n 

X\ e 1 e 1 > \ e 1 e 


Therefore, 


* Xi e e 


a.n j2nf.n a.n j2iTf.n 

1 - 1 > \ e 1 e 1 , where kX t = X' x 


Therefore, 

a', n jZitf.n a.n 32 itf,n a.n a.n 

X, e e > X, e 1 e , where ks = e 




(4 7 7) 


Therefore, 

a 'i n = a i n + ln k 

A g 

Note that X A are always positive, sinoe, they represent the power of the 
corresponding sinusoidal components Also, as cr e 4 and p are positive, it is evident 

ft 2 \ e 

from equation <4 7 6) that the term — p-*- will always be positive Clearly, the 
multiplying factor, k, in equation (4 7 7) will also be positive Moreover, it will 
always be grater than or equal to unity The later condition will be applicable for 
no noise conditions Thus, a\ will always be greater than a x , for any nonzero value 
of o e 2 

We observe that, the effect of noise is to increase the damping faotor It is 
well know that an increase of damping factor, in s domain, shifts the corresponding 
pole near the center of the unit circile, m z domain, l e , the pole m the z domain 
moves radially inward This radial inward shift is termed as radial bias From 
equation (4 7 6) we note that, in no noise oonditions, this radial shift is negligible 

Furthermore, we maintain that the angular displacement of the poles in z 
domain should theoretically be almost negligible This means that, idaaly, the 
frequency estimates must be very accurate, even for low SNRs 

Note that the entire left half of the s-plane is mapped in the unit circle in 
the z-plane Moreover, a pole at minus infinity m s-plane is equivalent to the pole 
at the oente of the unit circle in z-plane Thus, if the actual pole is away from the 
unit circle but inside it, then even a slight radial perturbation in its position in 
the z-plane will amount to a logarithmio change of its position in the s-plane Thus, 
while comparing simulation results, the position of the pole in the z-plane plays a 
very important role, since, the same percentage of change m the pole position, 
when it is away from the unit oirole, amounts to a very large error The unit 
circle repesentation, hence, is a somewhat elusive repesentation However, m 
applications, where the maximum and minimum damping factors differ widely, the unit 

•g. 

circle representation will be more convenient because of its concise nature 

Fig 4 14 shows two (of the four) signal pole estimates in the s domain for 
I the case depicted m Fig 4 6 Comparing Fig 4 14 and Fig 4 6 we observe that the 
I variations in the damping factor estimates are more apparent from Fig 4 14 than 
I from Fig 4 6 Thus, the s domain representation is more appropriate 



CHAPTER 5 

UNIFORM SAMPLING AS A SPECIAL CASE OF 
NONUNFQRM SAMPLING 


S I INTRODUCTION - 


In this chapter, we will consider the special case of parameter estimation 
from a uniformly sampled data The purpose of this chapter is to enlighten the 
noise cancelling features inherent m the orthogonal polynomial approximation 
procedure The effectiveness of the teohnique given in this thesis, will be clear 
from simulation results It is only expected that, same accuracy m pole estimates 
will be maintained for other similar cases 

5 2 EXAMPLE 1 - 


In this example, we will use equation (2 2 3) to approximate the given uniform 
samples This is to say that, Ct k k « 1, 2, , K3 will be uniformly distributed 

<Note that we are considering uniform sampling as a special case of nonuniform 
sampling) Thus, the white noise corrupted samples, y(t k ), will represent the given 
uniformly sampled corrupted signal, whose parameters are to be estimated In this 
example, we will consider the specific time-series, of r(t k >, given by, 

r(t k ) = e S±tk + A Z e S2tk + A 3 e S3tk + A 4 ***** + w(t k ) <52±> 


where, w(t k ) are zero mean independent, real, Gaussian stationary random variables, 

2 

with unknown variance, c w , and 


s L m 5 55 x 10 " 3 + j 0 OB, 

s 3 = 6 66 x iO " 3 + 3 0 li, 

t a 3 5, A 2 = 3 5, A 3 


S 2 - 5 55 x iQ- 3 - 3 0 08, 
s 4 = 6 66 x IO -3 - 3 0 11 , 
15, A 4 = 15, 


and, t k are uniformly spaced 

First, we preproc&ss y(t k ) using orthogonal polynomial approximation, to 
, generate the approximating signal, xCnl, as given by equation (4 2 2) Note that the 
SNR in xCnl is greater than the SNR in ytt^) Dnce the uniform samples, xCn3, n = 0, 
: 1, , L, are known, we use the estimator given by equation (4 6 4), to to find 







Eia.52 

letotrhOMfe Preprocessing. 

sm = £ da e = 2ql 

Uniform sampling. 





the estimate of the parameter vector, A MPC The SNR in y^) is 6 dB and tha moc * el 
order, p, used for simulation, is 20 Note that K * L * 50 The signal poles are 
calculated from the parameter vector, App^, by utilizing equation (4 2 3) The pole 
estimates for 20 different realizations of w(t k > are shown in Fig 5i 


5 3 EXAMPLE 2 - 


In this example, instead of approximating the given uniformly sampled data 
by y(t k ), we approxmate it by xtnD, as given by equation (4 2 2) The time series for 
rCnl, used in this example, was same as the one used in example i This approach, 
however, differs from the previous one m the sense that, we have not 
preprocessed the given data samples As xM is known we can directly use equation 
<4 6 4) to find the parameter vector, A MPc The SNR in xtnl is 6 dB and the model 
order, p, used for simulation, is 20 Also, L = 50 Note, however, that R ea« in ^ 1S 
case, will be an identity matrix, since, etnl, as given by equation (4 2 2), will 
represent a white noise sequence The pole estimates, m this case, - or 
different realizations of the white noise sequence, e£nl, are platted in Fig 5 2 


54 CONCLUSIONS - 


The table given below shows the estimates of the bias and variance of the 


signal parameters for the two oases mentioned above 





Comparing the bias and variance estimates for different signal parameters, 
Me can conclude that the bias and variance estimates in Example 1 are small as 
compared to those in Example 2 

Thus, preprocessing of data points is certainly helpful m increasing the 
accuracy of the estimates, although, the increase in the computaional burden 
should not be overlooked Lastly, we comment that, the estimation procedures used 
m Examples i and 2 gave similar performance, above a SNR level qf 15 dB Thus, 
for accurate parameter estimates, the additional complexity of preprocesing will 
be advantageous when the SNR level is below this threshold of 15 dB 
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APPENDIX A 


From equation <2 3 9), we can write (P^Pf* as, 


<P T P) _1 « 


t . *>o Z V 


k a i 


k «* i 


0l\> 
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k ■ 1 


PN-l^V 


(Ai) 


From equation <2 3 6), P T is given bw, 
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ppdjL) 
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PpCt^) 
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p i <t 2 ) 

D L (t K ) 
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Using equation Ai and A2, <P^Pf * P^*, pan be expressed as, 
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Using equations (A3) and (3 2 5), we can write q n m» which is the element at the nth 

T T 

row and mth column of the matrix P'<P P) P , as, 

Pnrn = | — | , where p^tnl = p i _j,CnT3 <A4) 

1 - 1 S Pi-ltt^l 

k = 1 

Note that the transformation coefficients, q n m» will be correlated with each other 
This comes from the fact that the approximating polynomial is a deterministic 
signal Thus, the value of this polynomial at any time instant tn will be correlated 
with the value of the same polynomial at different time instant, say t n+ j < It is 
apperent, that the amount of correlation will decrease as the lag, k, is increased 


APPENDIX B 


Let xCn3, be a noise corrupted signal, such that, it consists of the signal, sCnl, and 
the noise, w[n3, respectively, 1 e , 

xtn3 = stn3 + wtn3 


Let xC03, , xtp3, xCp+13, , xCL-13, be the L uniformly spaced samples of the 

noise corrupted process xtn3 Thus, we can write the following equations, 


xCp3 

a s[p3 

+ etp3 

xCp+13 

= SCp+13 

+ eCp+13 

xCL-13 

= sCL-i] 

+ aCL-13 


Equation (Bi) can be written in matrix notation as, 


X = 

S + 

W 


B2 

where, 

X = 

C xtp] 

xCp+13 

xCL-13 ] T 

B3 

S = 

C slp3 

sCp+13 

sCL-13 ] T 

B4 

W a 

[ wtp3 

wtp+13 

wCL-13 ] T 

B5 


Let the autocorrelation matrix, Rxx* he defined as 
R xx = E(XX T ) 


Therefore, 


Rxx = 


x[p3 

xCp+13 


xCL-13 



xCp+13 


xCL-iD^J 


BE 


E 


Therefore, 

£<x[p3xCp3> 

£(x[p3x£p+13) 

E(xCp3xCL-J,3) 



E(xCp+i3xCp3) 

£<x[p+ilxCp+13> 

£(xCp+i3xCL-13) 


Rxx “ 

£<xCL-JL]x[p3> 

E<xCL-i3x[p+i3) 

E(xCL-13xCL-13) 

B7 


where, the autocorrelation function, rxxi is given by, 

r xx [ n+k, rt * E( x*[n+k3 xtn3 j B8 

« £ ( ( s*Cn+kD + w*Cn+kl ) ( stn] + wCn3 ) ) 

As the orosscorrelation functions, 


rsw^n+k, n3 » £ ( s*Cn+kD wtn3 ) 

B9 

r WB [n+k, nD = E ( w*Cn+k3 stn3 ] 

BIO 

zero, we get, 


rxxtn+k, nl = E ( s*Cn+k] sCn3 ] + E ( w*[n+k] wCn3 j 


= r ss [n+k, n3 + r ww [n+k, nl 

Bli 


Substituting, equation Bil in B7, we get, 


R xx = Rss + R ww 


B12 


If sCnl and wCn3 are WSS, then Rxx is a symmetric Toeplitx matrix If s[nl and 
wCnl are looally WSS [43, then Rxx is a symmetnx matrix, but not Toeplitz If the 
processes are nonstationary C43 then, Rxx has np speoial properties, such as 
symmetricity or Toeplitz structure 


APPENDIX C 


Let xCnl be a process, which is assumed to be locally WSS C43 Let xtQ3, , xCL-13, 
be L uniformly sampled samples of xCnl We define a matrix X, such that, 


X 


xCp-13 

xCpl 


<L-p) 


1/2 


xCp-23 

xCp-13 


xtL-23 xCL-33 


xCQ3 

xElD 


xCL-p-13 


Thus, the p by p autocorrelation matrix, R xx > can be written as, 
R xx = X H X 


Cl 


C2 


Let the singular value decomposition (SVD) of X be given by, 


X = U V 


iH 


where, 

U = [ u A u 2 


u L-p ^ 


V = [ Vj_ v 2 


and 


£ = 


0 


0 

0 

0 

0 

0 


X, 0 


0 

0 


0 0 


Vp ] 


0 


0 

0 


D 

X p 

0 


C3 

C4 

C5 


C6 


i 




c ub_t-i* jU.g eq<i3**on tCw' in equation <CD, we get, 


p,< - i U Z v H ) H ( u 2 v H I 
= V £ T iU H U) £ V H 

= v v H r vr v H 

where 



From equation <C7) we observe that the column and row eigen vectors are same 


i 
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